home *** CD-ROM | disk | FTP | other *** search
/ Amiga Plus Leser 15 / Amiga Plus Leser CD 15.iso / Tools / Development / yacas_alg / yacas_morphos / share / yacas / standard.ys < prev    next >
Encoding:
Text File  |  2002-03-13  |  12.9 KB  |  546 lines
  1.  
  2. /* See the documentation on the assignment of the precedence of the rules.
  3.  */
  4.  
  5. /* Some very basic functions that are used always any way... */
  6.  
  7.  
  8. /* Implementation of Nth that allows extending. */
  9. RuleBase("Nth",{alist,aindex});
  10. Rule("Nth",2,10,
  11.     MathAnd(Equals(IsFunction(alist),True),
  12.             Equals(IsInteger(aindex),True),
  13.             Not(Equals(Head(Listify(alist)),Nth))
  14.             ))
  15.      MathNth(alist,aindex);
  16.  
  17.  
  18.  
  19.  
  20. Rule("Nth",2,14,
  21.      MathAnd(Equals(IsString(alist),True),IsList(aindex))
  22.     )
  23. [
  24.   Local(result);
  25.   result:="";
  26.   ForEach(i,aindex) [ result := result : StringMid(i,1,alist); ];
  27.   result;
  28. ];
  29.  
  30. Rule("Nth",2,15,Equals(IsString(alist),True))
  31. [
  32.   StringMid(aindex,1,alist);
  33. ];
  34.  
  35.  
  36. Rule("Nth",2,20,Equals(IsList(aindex),True))
  37. [
  38.   Map({{ii},alist[ii]},{aindex});
  39. ];
  40.  
  41. Rule("Nth",2,30,
  42.    MathAnd(
  43.            Equals(IsGeneric(alist),True),
  44.            Equals(GenericTypeName(alist),"Array"),
  45.            Equals(IsInteger(aindex),True)
  46.           )
  47.     )
  48. [
  49.   ArrayGet(alist,aindex);
  50. ];
  51.  
  52.  
  53.  
  54. Rule("Nth",2,40,Equals(IsString(aindex),True))
  55. [
  56.   Local(as);
  57.   as := Assoc(aindex,alist);
  58.   If (Not(Equals(as,Empty)),Set(as,Nth(as,2)));
  59.   as;
  60. ];
  61.  
  62.  
  63.  
  64. Function("NrArgs",{aLeft}) Length(Listify(aLeft))-1;
  65.  
  66. 10 # IsNonObject(Object(_x)) <-- False;
  67. 20 # IsNonObject(_x)         <-- True;
  68.  
  69. 1 # Numer(_x / _y)      <-- x;
  70. 2 # Numer(x_IsNumber)   <-- x;
  71. 1 # Denom(_x / _y)      <-- y;
  72. 2 # Denom(x_IsNumber)   <-- 1;
  73.  
  74. /* Cache the most recently computed value of Pi and its precision
  75.   This helps because right now every N() requires a new evaluation of Pi
  76.  */
  77. Set(PiCache, {40, 3.1415926535897932384626433832795028841971694});
  78.     /* return numerical value recalculated only if necessary
  79.       we could get more speed if we also passed the previous value of pi
  80.       to MathPi() as its initial guess but that can wait for now */
  81. Function("Pi", {})
  82. [
  83.     Local(new'prec, new'pi);
  84.     Set(new'prec, GetPrecision());
  85.     If ( LessThan(MathNth(PiCache,1), new'prec),
  86.     [
  87.         Set(new'pi, MathPi());
  88.         Set(PiCache, {new'prec, new'pi});
  89.         new'pi;
  90.     ],
  91.     MathNth(PiCache,2)
  92.     );
  93. ];
  94.  
  95. Set(Numeric,False);
  96. Clear(Pi);
  97. Function("N",{aNumberBody})
  98. [
  99.   Local(prevNumeric,prevPi,numericresult);
  100.   Set(prevNumeric,Numeric);
  101.   Set(prevPi,Pi);
  102.   Set(Pi,Pi());
  103.   Set(Numeric,True);
  104.   Set(numericresult,Eval(aNumberBody));
  105.   Set(Numeric,prevNumeric);
  106.   Set(Pi,prevPi);
  107.   numericresult;
  108. ];
  109. UnFence("N",1);
  110.  
  111. Function("N",{aNumberBody,aDigits})
  112. [
  113.   Local(prevNumeric,prevPi,numericresult);
  114.   Set(prevNumeric,Numeric);
  115.   Set(prevPi,Pi);
  116.   Set(prevdigits,GetPrecision());
  117.   Precision(aDigits);
  118.   Set(Pi,Pi());
  119.   Set(Numeric,True);
  120.   Set(numericresult,Eval(aNumberBody));
  121.   Set(Numeric,prevNumeric);
  122.   Set(Pi,prevPi);
  123.   Precision(prevdigits);
  124.   numericresult;
  125. ];
  126. UnFence("N",2);
  127.  
  128.  
  129. Set(Verbose,False);
  130. Function("V",{aNumberBody})
  131. [
  132.   Local(prevVerbose,result);
  133.   Set(prevVerbose,Verbose);
  134.   Set(Verbose,True);
  135.   Set(result,Eval(aNumberBody));
  136.   Set(Verbose,prevVerbose);
  137.   result;
  138. ];
  139. HoldArg("V",aNumberBody); 
  140. UnFence("V",1);
  141.  
  142. Function("++",{aVar})
  143. [
  144.    MacroSet(aVar,MathAdd(Eval(aVar),1));
  145. ];
  146. UnFence("++",1);
  147. HoldArg("++",aVar);
  148.  
  149.  
  150. Function("--",{aVar})
  151. [
  152.    MacroSet(aVar,MathSubtract(Eval(aVar),1));
  153. ];
  154. UnFence("--",1);
  155. HoldArg("--",aVar);
  156.  
  157.  
  158. Function("TableForm",{list})
  159. [
  160.   Local(i);
  161.   ForEach(i,list)
  162.   [
  163.     Write(i);
  164.     NewLine();
  165.   ];
  166.   True;
  167. ];
  168.  
  169.  
  170.  
  171. /* Standard arithmetic */
  172.  
  173.  
  174. /* Addition */
  175.  
  176. 100 # + _x  <-- x;
  177.  
  178. 50 # x_IsNumber + y_IsNumber <-- MathAdd(x,y);
  179.  
  180. 100 # 0 + _x    <-- x;
  181. 100 # _x + 0    <-- x;
  182. 100 # _x + _x   <-- 2*x;
  183. 101 # _x + - _y <-- x-y;
  184. 101 # _x + (- _y)/(_z) <-- x-(y/z);
  185. 101 # (- _y)/(_z) + _x  <-- x-(y/z);
  186. 101 # (- _x) + _y <-- y-x;
  187. 102 # _x + y_IsNegativeNumber <-- x-(-y);
  188. 102 # _x + y_IsNegativeNumber * _z <-- x-((-y)*z);
  189. 102 # _x + (y_IsNegativeNumber)/(_z) <-- x-((-y)/z);
  190. 102 # (y_IsNegativeNumber)/(_z) + _x  <-- x-((-y)/z);
  191. 102 # (x_IsNegativeNumber) + _y <-- y-(-x);
  192. 150 # _n1 / _d + _n2 / _d <-- (n1+n2)/d;
  193.  
  194. 200 # (x_IsNumber + _y)_Not(IsNumber(y)) <-- y+x;
  195. 200 # ((_y + x_IsNumber) + _z)_Not(IsNumber(y) Or IsNumber(z)) <-- (y+z)+x;
  196. 200 # ((x_IsNumber + _y) + z_IsNumber)_Not(IsNumber(y)) <-- y+(x+z);
  197. 200 # ((_x + y_IsNumber) + z_IsNumber)_Not(IsNumber(x)) <-- x+(y+z);
  198.  
  199. 210 # x_IsNumber + (y_IsNumber / z_IsNumber) <--(x*z+y)/z;
  200. 210 # (y_IsNumber / z_IsNumber) + x_IsNumber <--(x*z+y)/z;
  201. 210 # (x_IsNumber / v_IsNumber) + (y_IsNumber / z_IsNumber) <--(x*z+y*v)/(v*z);
  202.  
  203.  
  204. LocalSymbols(x)
  205. [
  206.   220 # + x_IsList          <-- MapSingle("+",x);
  207. ];
  208. 220 # xlist_IsList + ylist_IsList <-- Map("+",{xlist,ylist});
  209.  
  210. SumListSide(_x, y_IsList) <--
  211. [
  212.    Local(i,result);
  213.    result:={};
  214.    For(i:=1,i<=Length(y),i++)
  215.    [ DestructiveInsert(result,i,x + y[i]); ];
  216.    result;
  217. ];
  218.  
  219. 240 # (x_IsList + _y)_Not(IsList(y)) <-- SumListSide(y,x);
  220. 241 # (_x + y_IsList)_Not(IsList(x)) <-- SumListSide(x,y);
  221.  
  222. 250 # Infinity + y_IsNumber <-- Infinity;
  223. 250 # x_IsNumber + Infinity <-- Infinity;
  224.  
  225.  
  226. /* Subtraction arity 1 */
  227.  
  228. //50 # -0 <-- 0;
  229. 51 # -Undefined <-- Undefined;
  230. 54 # - (- _x)      <-- x;
  231. 55 # (- (x_IsNumber)) <-- MathSubtract(x);
  232. 110 # - (_x - _y) <-- y-x;
  233. 111 # - (x_IsNumber / _y) <-- (-x)/y;
  234. LocalSymbols(x)
  235. [
  236.   200 # - (x_IsList) <-- MapSingle("-",x);
  237. ];
  238.  
  239. /* Subtraction arity 2 */
  240. 50  # x_IsNumber - y_IsNumber <-- MathSubtract(x,y);
  241. 50  # x_IsNumber - y_IsNumber <-- MathSubtract(x,y);
  242. 60  # Infinity - Infinity <-- Undefined;
  243. 100 # 0 - _x <-- -x;
  244. 100 # _x - 0 <-- x;
  245. 100 # _x - _x <-- 0;
  246.  
  247. 110 # _x - (- _y) <-- x + y;
  248. 110 # _x - (y_IsNegativeNumber) <-- x + (-y);
  249. 111 # (_x + _y)- _x <-- y;
  250. 111 # (_x + _y)- _y <-- x;
  251. 112 # _x - (_x + _y) <-- - y;
  252. 112 # _y - (_x + _y) <-- - x;
  253. 113 # (- _x) - _y <-- -(x+y);
  254. 113 # (x_IsNegativeNumber) - _y <-- -((-x)+y);
  255. 113 # (x_IsNegativeNumber)/_y - _z <-- -((-x)/y+z);
  256.  
  257. /* TODO move to this precedence everywhere? */
  258. 10 # (x_IsList) - (y_IsList) <-- Map({{xarg,yarg},xarg-yarg},{x,y});
  259.  
  260. 240 # (x_IsList - y_IsNonObject)_Not(IsList(y)) <-- -(y-x);
  261.  
  262. 241 # (x_IsNonObject - y_IsList)_Not(IsList(x)) <--
  263. [
  264.    Local(i,result);
  265.    result:={};
  266.    For(i:=1,i<=Length(y),i++)
  267.    [ DestructiveInsert(result,i,x - y[i]); ];
  268.    result;
  269. ];
  270.  
  271. 250  # Infinity - y_IsNumber <-- Infinity;
  272. 250  # x_IsNumber - Infinity <-- - Infinity;
  273.  
  274.  
  275. 210 # x_IsNumber - (y_IsNumber / z_IsNumber) <--(x*z-y)/z;
  276. 210 # (y_IsNumber / z_IsNumber) - x_IsNumber <--(y-x*z)/z;
  277. 210 # (x_IsNumber / v_IsNumber) - (y_IsNumber / z_IsNumber) <--(x*z-y*v)/(v*z);
  278.  
  279.  
  280. /* Multiplication */
  281.  
  282. 50  # x_IsNumber * y_IsNumber <-- MathMultiply(x,y);
  283. 100 #  1  * _x  <-- x;
  284. 100 # _x  *  1  <-- x;
  285. 100 # (_f  * _x)_(f= -1)  <-- -x;
  286. 100 # (_x  * _f)_(f= -1)  <-- -x;
  287.  
  288. 100 # x_IsNumber * (y_IsNumber * _z) <-- (x*y)*z;
  289. 100 # x_IsNumber * (_y * z_IsNumber) <-- (x*z)*y;
  290.  
  291. 100 # (_x * _y) * _y <-- x * y^2;
  292. 100 # (_x * _y) * _x <-- y * x^2;
  293. 100 # _y * (_x * _y) <-- x * y^2;
  294. 100 # _x * (_x * _y) <-- y * x^2;
  295. 100 # _x * (_y / _z) <-- (x*y)/z;
  296. 100 # (_y / _z) * _x <-- (x*y)/z;
  297. 100 # (_x * y_IsNumber)_Not(IsNumber(x)) <-- y*x;
  298.  
  299. 100 # (_x * _y)* _x ^ n_IsInteger <-- y * x^(n+1);
  300. 100 # (_y * _x)* _x ^ n_IsInteger <-- y * x^(n+1);
  301.  
  302. 105 # x_IsNumber * -(_y) <-- (-x)*y;
  303. 105 # (-(_x)) * (y_IsNumber) <-- (-y)*x;
  304.  
  305. 106 # _x * -(_y) <-- -(x*y);
  306. 106 # (- _x) * _y <-- -(x*y);
  307.  
  308. 107 # -( (-(_x))/(_y)) <-- x/y;
  309. 107 # -( (_x)/(-(_y))) <-- x/y;
  310.  
  311.  
  312.  
  313. 200 # x_IsMatrix * y_IsMatrix <-- 
  314. [
  315.    Local(i,j,k,row,result);
  316.    result:=ZeroMatrix(Length(x),Length(y[1]));
  317.    For(i:=1,i<=Length(x),i++)
  318.    For(j:=1,j<=Length(y),j++)
  319.    For(k:=1,k<=Length(y[1]),k++)
  320.    [
  321.      row:=result[i];
  322.      row[k]:= row[k]+x[i][j]*y[j][k];
  323.    ];
  324.    result;
  325. ];
  326.  
  327.  
  328. 210 # x_IsMatrix * y_IsList <-- 
  329. [
  330.    Local(i,result);
  331.    result:={};
  332.    For(i:=1,i<=Length(x),i++)
  333.    [ DestructiveInsert(result,i,x[i] . y); ];
  334.    result;
  335. ];
  336.  
  337.  
  338. 240 # (x_IsList * y_IsNonObject)_Not(IsList(y)) <-- y*x;
  339. 241 # (x_IsNonObject * y_IsList)_Not(IsList(x)) <--
  340. [
  341.    Local(i,result);
  342.    result:={};
  343.    For(i:=1,i<=Length(y),i++)
  344.    [ DestructiveInsert(result,i,x * y[i]); ];
  345.    result;
  346. ];
  347.  
  348. 249  # 0 * Infinity <-- Undefined;
  349. 250  # x_IsNumber * y_IsInfinity <-- Sign(x)*y;
  350. 250  # x_IsInfinity * y_IsNumber <-- Sign(y)*x;
  351.  
  352.  
  353. /* Note: this rule MUST be past all the transformations on
  354.  * matrices, since they are lists also.
  355.  */
  356. 230 # aLeft_IsList * aRight_IsList <-- Map("*",{aLeft,aRight});
  357. 242 # (x_IsInteger / y_IsInteger) * (v_IsInteger / w_IsInteger) <-- (x*v)/(y*w);
  358. 243 #  x_IsInteger * (y_IsInteger / z_IsInteger) <--  (x*y)/z;
  359. 243 #  (y_IsInteger / z_IsInteger) * x_IsInteger <--  (x*y)/z;
  360.  
  361. 400 # _x * _x <-- x^2;
  362. 401 # 0 * _x <-- 0;
  363. 401 # _x * 0 <-- 0;
  364.  
  365.  
  366. /* Division */
  367.  
  368. 50 # _x / 0 <-- Undefined;
  369. 51 # 0 / _x <-- 0;
  370. 52 # (x_IsNumber / y_IsNumber)_((Numeric = True) Or
  371.                                 Not(IsInteger(x) And IsInteger(y))) <--
  372.     MathDivide(x,y);
  373.  
  374.  
  375. 51 # x_IsNumber / y_IsNegativeNumber <-- (-x)/(-y);
  376.  
  377. 51 # (x_IsNonZeroInteger / y_IsNonZeroInteger)_(MathGcd(x,y) > 1) <--
  378.      [
  379.        Local(gcd);
  380.        Set(gcd,MathGcd(x,y));
  381.        MathDivide(x,gcd)/MathDivide(y,gcd);
  382.      ];
  383.     
  384.  
  385. 90 # x_IsInfinity / y_Infinity <-- Undefined;
  386. 91  # x_IsInfinity / y_IsNumber <-- Sign(y)*x;
  387.  
  388.  
  389. 100 # _x / _x <-- 1;
  390. 100 # _x /  1 <-- x;
  391. 100 # (_x / y_IsNegativeNumber) <-- -x/(-y);
  392. 100 # (_x / - _y) <-- -x/y;
  393. 200 # (_x / _y)/ _z <-- x/(y*z);
  394.  
  395. 230 # _x / (_y / _z) <-- (x*z)/y;
  396. 240 # xlist_IsList / ylist_IsList <-- Map("/",{xlist,ylist});
  397.  
  398.  
  399. 250 # x_IsList / _y <--
  400. [
  401.    Local(i,result);
  402.    result:={};
  403.    For(i:=1,i<=Length(x),i++)
  404.    [ DestructiveInsert(result,i,x[i] / y); ];
  405.    result;
  406. ];
  407.  
  408. 250 # _x / y_IsList <--
  409. [
  410.    Local(i,result);
  411.    result:={};
  412.    For(i:=1,i<=Length(y),i++)
  413.    [ DestructiveInsert(result,i,x/y[i]); ];
  414.    result;
  415. ];
  416.  
  417. 250 # _x / Infinity <-- 0;
  418. 250 # _x / (-Infinity) <-- 0;
  419.  
  420.  
  421. 400 # 0 / _x <-- 0;
  422.  
  423. /* Faster version of raising power to 0.5 */
  424. 50 # (x_IsPositiveNumber ^ (1/2))_(Numeric) <-- Sqrt(x);
  425. 50 # (x_IsPositiveNumber ^ (1/2))_IsInteger(MathSqrt(x)) <-- MathSqrt(x);
  426. 58 # 1 ^ n_IsInfinity <-- Undefined; 
  427. 59 # _x ^ 1 <-- x;
  428. 59 # 1 ^ _n <-- 1;
  429. 59 # 0 ^ 0 <-- Undefined;
  430. 60 # 0 ^ n_IsRationalOrNumber <-- 0;
  431.  
  432. /* Regular raising to the power. */
  433. 61 # Infinity ^ (y_IsNegativeNumber) <-- 0;
  434. 61 # (-Infinity) ^ (y_IsNegativeNumber) <-- 0;
  435. 61 # x_IsPositiveNumber ^ y_IsPositiveNumber <-- MathPower(x,y);
  436. 61 # x_IsPositiveNumber ^ y_IsNegativeNumber <-- (1/MathPower(x,-y));
  437.  
  438. 90 # (-_x)^m_IsEven <-- x^m;
  439. 90 # (Sqrt(x_IsConstant))_(IsNegativeNumber(N(x))) <-- Complex(0,Sqrt(-x));
  440. 91 # (x_IsConstant ^ (m_IsOdd / p_IsOdd))_(IsNegativeNumber(Re(N(x)))) <--
  441.      -((-x)^(m/p));
  442. 92 # (x_IsNegativeNumber ^ y_IsNumber)_Numeric <-- Exp(y*Ln(x));
  443.  
  444.  
  445. 70  # (_x ^ m_IsRationalOrNumber) ^ n_IsRationalOrNumber <-- x^(n*m);
  446.  
  447. 80 # (x_IsNumber/y_IsNumber) ^ n_IsPositiveInteger <-- x^n/y^n;
  448. 80 # x_IsNegativeNumber ^ n_IsEven <-- (-x)^n;
  449. 80 # x_IsNegativeNumber ^ n_IsOdd <-- -((-x)^n);
  450.  
  451. /*TODO remove?
  452. 80 # x_IsRationalOrNumber ^ n_IsPositiveInteger <--
  453. [
  454.   Local(mult,result);
  455.   Set(mult,x);
  456.   Set(result,1);
  457.   While(n>0)
  458.   [
  459.     If(n&1 != 0, Set(result, result*mult));
  460.     Set(n,n>>1);
  461.     Set(mult,mult*mult);
  462.   ];
  463.   result;
  464. //  x*x^(n-1);
  465. ];
  466. */
  467.  
  468. 100  # ((_x)*(_x ^ _m)) <-- x^(m+1);
  469. 100  # ((_x ^ _n)*(_x ^ _m)) <-- x^(m+n);
  470. 100  # ((_x)^(_n))*((_y)^(_n)) <-- (x*y)^n;
  471.  
  472. 100  # ((x_IsNumber)^(n_IsNumber/(_m)))_(n>1) <-- MathPower(x,n)^(1/m);
  473.  
  474. 100 # Sqrt(_n)^(m_IsEven) <-- n^(m/2);
  475.  
  476.  
  477. /*
  478. 100 # ((_x)^(_n))^((_m)/(_p)) <-- x^(n*m/p);
  479. 100 # ((_x)^((_m)/(_p)))^(_n) <-- x^(n*m/p);
  480.  
  481. 100 # ((_x)^((_m)/(_p)))^((_n)/(_q)) <-- x^((n*m)/(p*q));
  482.  */
  483. 200 # xlist_IsList ^ nlist_IsList <-- Map("^",{xlist,nlist});
  484. 201 # xlist_IsList ^ n_IsConstant <-- Map({{xx},xx^n},{xlist});
  485. 202 # _x ^ n_IsList <-- Map({{xx},x^xx},{n});
  486.  
  487. 249 # x_IsInfinity ^ 0 <-- Undefined;
  488. 250 # Infinity ^ (_n) <-- Infinity;
  489. 250 # ((-Infinity) ^ (n_IsNumber))_(IsEven(n)) <-- Infinity;
  490. 250 # ((-Infinity) ^ (n_IsNumber))_(IsOdd(n)) <-- -Infinity;
  491.  
  492. 250 # (x_IsNumber ^ Infinity)_(x> -1 And x < 1) <-- 0;
  493. 250 # (x_IsNumber ^ Infinity)_(x< -1) <-- - Infinity;
  494. 250 # (x_IsNumber ^ Infinity)_(x> 1) <-- Infinity;
  495.  
  496. 250 # (x_IsNumber ^ -Infinity)_(x> -1 And x < 1) <-- Infinity;
  497. 250 # (x_IsNumber ^ -Infinity)_(x< -1) <-- 0;
  498. 250 # (x_IsNumber ^ -Infinity)_(x> 1) <-- 0;
  499.  
  500.  
  501. 400 # _x ^ 0 <-- 1;
  502.  
  503.  
  504. RuleBase("==",{left,right});
  505. RuleBase("!==",{left,right});
  506.  
  507.  
  508. Function("+-",{x,y})  x + -y;
  509. Function("/-",{x,y})  x / -y;
  510. Function("*-",{x,y})  x * -y;
  511. Function("^-",{x,y})  x ^ -y;
  512.  
  513. Function(":=-",{xleft,yminright}) Apply(":=",{xleft,-yminright}); 
  514. HoldArg(":=-",xleft);
  515. HoldArg(":=-",yminright);
  516.  
  517. Function(":=+",{xleft,yright}) Apply(":=",{xleft,yright}); 
  518. HoldArg(":=+",xleft);
  519. HoldArg(":=+",yright);
  520.  
  521.  
  522. a_IsPositiveInteger & b_IsPositiveInteger <-- BitAnd(a,b);
  523. a_IsPositiveInteger | b_IsPositiveInteger <-- BitOr(a,b);
  524. a_IsPositiveInteger % b_IsPositiveInteger <-- Mod(a,b);
  525.  
  526. RuleBase("if",{predicate,body});
  527. (if(True) _body) <-- Eval(body);
  528. HoldArg("if",body);
  529. UnFence("if",2);
  530.  
  531. RuleBase("else",{ifthen,otherwise});
  532. 0 # (if (_predicate) _body else _otherwise)_(Eval(predicate) = True) <--
  533.      Eval(body);
  534. 0 # (if (_predicate) _body else _otherwise)_(Eval(predicate) = False) <--
  535.      Eval(otherwise);
  536. 1 # (if (_predicate) _body else _otherwise) <--
  537.     UnList({Atom("else"),
  538.             UnList({Atom("if"), (Eval(predicate)), body}),
  539.             otherwise});
  540. HoldArg("else",ifthen);
  541. HoldArg("else",otherwise);
  542. UnFence("else",2);
  543.  
  544.  
  545.  
  546.